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ABSTRACT 

We present a study of the cosmological Lya emission signal at z > 4. Our goal is to 
predict the power spectrum of the spatial fluctuations that could be observed by an 
intensity mapping survey. The model uses the latest data from the HST legacy fields 
and the abundance matching technique to associate UV emission and dust proper¬ 
ties with the halos, computing the emission from the interstellar medium (ISM) of 
galaxies and the intergalactic medium (IGM), including the effects of reionization, 
self-consistently. The Lya intensity from the diffuse IGM emission is 1.3 (2.0) times 
more intense than the ISM emission at z = 4(7); both components are fair tracers of 
the star-forming galaxy distribution. However the power spectrum is dominated by 
ISM emission on small scales (k > 0.01/iMpc -1 ) with shot noise being significant only 
above k = 1/iMpcG 1 . At very lange scales (k < 0.01/iMpc -1 ) diffuse IGM emission 
becomes important. The comoving Lya luminosity density from IGM and galaxies, 
P l Lya — 8.73(6.51) x 10 4o erg s _1 Mpc -3 and = 6.62(3.21) x 10 4o erg s -1 Mpc -3 
at z = 4(7), is consistent with recent SDSS determinations. We predict a power 
k 3 P Lya (k,z)/2TT 2 = 9.76 x 10 _4 (2.09 x l(U 5 )nW 2 m- 4 sr- 2 at 2 = 4(7) for k = 
0.1/iMpc -1 . 

Key words: cosmology: observations - intergalactic and interstellar medium - inten¬ 
sity mapping - large-scale structure of universe 


1 INTRODUCTION 

According to the most popular cosmological framework, the 
ACDM model, the astonishing diversity present in the local 
universe originated from tiny density fluctuations in a ho¬ 
mogeneous hot plasma. During the Hubble expansion, the 
dark matter (DM) in the most overdense regions started to 
collapse through gravitational instabilities, forming virial- 
ized halos where baryons could cool and form stars. This 
process gave birth to primordial galaxies and, with billions 
of years, every structure we can observe today. The impact 
of these early objects on the subsequent cosmic evolution, 
mediated by a number of feedback processes, has been dra¬ 
matic; in addition, their radiative energy input powered the 
last phase transition in the universe during the Epoch of 
Reionization (EoR, (Barkana fe Loebp OOl l). 

It is quite surprising that, for fifty years, we have been 
able to observe directly both the beginning, through the 
CMB, and the ending of this process. However the study of 
the ancient universe proved to be extremely challenging and 
only in the last decade we have been able to explore the EoR 
latest stages. 

One of the main challenges is that the first galaxies are 
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very faint and hard to detect even with our most power¬ 
ful telescopes. Up to date, the deepest “drop-out” photo¬ 
metric surveys carried out by the Hubble Space Telescope 
(HST), could detect ~ 700 galaxies at z > 8 (Bouwens et al. 
|2014a|. Moreover, such sources are the most massive outliers 
of the galaxy population and therefore bear little informa¬ 
tion about the sources producing the bulk of the ionizing 
radiation (Salvaterra et al. [2011). This difficulty will likely 
not be overcome by even the next generation of deep sur¬ 
veys, such as the JWST one. 

For this reason, a different strategy must be designed, 
which implies to forego the detection of individual sources 
and probe directly their large scale distribution. This idea 
can be implemented through the Intensity Mapping (IM) of 
selected emission lines (|Visbal fe~ Locb 2010. Visbal~et al.| 


2011). Each point in space is identified by an angular coor¬ 


dinate and a redshift, that can be measured knowing the fre¬ 
quency of the emitted and detected photons. If the relevant 
foregrounds can be removed, one can measure the distribu¬ 
tion of the cumulative galaxy emission in coarse voxels. For 
example, if the galaxy-galaxy correlation length is 1 Mpc 
and we are interested in structures larger than 10 Mpc, we 
can probe the sky in 10 Mpc voxels, corresponding to a spec¬ 
tral (angular) resolution \/8\ ss 300 (6 « 4') at z = 7, by 
grouping together the signal of ~ 1000 galaxies. 
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Several emission lines have been proposed for this kind 


of survey: 21cm (Furlanetto et al. 2006), CO 

Lidz et al. 

2011 

Righi et al. 2008 Breysse et al. 2014), H 2 ( 

Gong et al. 

2013 

I, ICIIJ 

1, Hell ( 

Gong et al.|2012 

Silva et al. |||2014 

Yue et al. 

2015 

Visbal et al.|2015 

1 and, of course, 

Lya (Silva 

et al.| 2013( |Pullen et al. 20141. Noticeably, each line probes 


different physical processes and, therefore, bears comple¬ 
mentary information on the sources. 

The Lya line most likely represents the optimal feature 
for a IM experiment. Historically, it has constantly been used 
for high-z galaxy searches as it is the most luminous UV line. 
Observers have undertaken deep surveys to detect high red- 
shift Lya emitters (LAE) ( |Ouchi et al.|20d8||2010||Matthee| 
et al.| [2015). These searches are affected by the same prob¬ 
lems mentioned above in the context of drop-out searches, 
with the additional complication arising from the fact that 
intergalactic H i can scatter the bulk of Lya photons out of 
the line of sight, making systematic detections of LAE dur¬ 


ing the EoR very challenging (Dayal et al. 2011). IM can 


overcome these problems, thanks to its sensitivity to even 
diffuse emission from the IGM. Therefore Lya IM seems a 
very promising tool to study the properties of early, faint 
and distributed EoR sources. 


Silva et al. (20131 studied this problem with seminu- 


merical tools, focusing on the EoR emission at z = 7. In a 


following work, Gong et al. (2014) tackled the problem of 


line confusion, proposing the masking of contaminated pix¬ 
els as a cleaning technique. In these works it is found that 
recombinations in the interstellar medium (ISM) of galax¬ 
ies largely dominate the Lya intensity and power spectrum 
(PS). Finally, Pullen et al.1 (2014] developed a simple analyt¬ 
ical model to study the evolution of the Lya power spectrum 
at z > 2: their results are qualitatively different from |SiIva| 
et al. (2013), in that they conclude that diffuse IGM emis¬ 


sion is the dominant source. Recently Croft et al. (2015) 


attempted to observe the large scale clustering of Lya emis¬ 
sion at z = 2 — 3 by cross-correlating the residuals in the 
SDSS spectra with QSOs. They claim a detection of a mean 
Lya surface brightness Xi 10 times more intense than the 
one inferred from LAE surveys, but still compatible with 
the unobscured Lya emission expected from LBGs. 

In this work we develop a detailed analytical model for 
the Lya emission and PS. We will use the latest data from 
the high redshift LBG in the HST legacy fields to compute 
the emission of ionizing and UV radiation and dust obscu¬ 
ration. Then we will model the physics of the ISM and of 
the IGM to obtain a data-constrained reionization history 
and expected Lya intensity. Our main aim is to extract the 
physical information encoded in a Lya intensity map, and 
assess whether it can be used to probe the earliest genera¬ 
tion of galaxies, their ISM and also the structure of the IGM 
during the EoR. 

The paper is organized as follows: in S|2]we introduce our 
model for star formation and in m we compute the reion¬ 
ization history accordingly; in SI we review the physics in¬ 
volved in Lya emission and we use it to predict the mean 
Lya intensity from the different processes involved. Finally, 
in (J5]we derive the power spectrum of the Lya fluctuations 
and discuss the results. We assume a flat ACDM cosmol¬ 
ogy compatible with the latest Planck results (h = 0.677, 
Qm = 0.31, n b = 0.049, U A = 1 - f2m, n = 0.97, a s = 0.82, 


2 HIGH REDSHIFT GALAXIES 

The first step towards the computation of the Lya inten¬ 
sity fluctuations is to accurately model the UV emissivity of 
high-redshift galaxies and their contribution to IGM reion¬ 
ization. This is because Lya emission mostly arises as a re¬ 
sult of recombinations in gas ionized by young, massive stars. 
Ionizing photons are produced also by Active Galactic Nu¬ 
clei; however, as there is a general consensus that reioniza¬ 


tion is largely driven by stars (but see Giallongo et al. (20151 


and Weigel et al. (2015 l), we will neglect in the following the 


presence of AGNs. In the following we model the UV emis¬ 
sivity of galaxies using an empirical approach based on the 
most recent UV Luminosity Functions measurements. 

2.1 UV emissivity 

As Lya emission is intimately related to massive stars which 
have a lifetime (few Myr) short compared to the Hubble 
time, it effectively depends only on the instantaneous Star 
Formation Rate (SFR), ip, of the galaxy. We seek a relation 
between the halo mass, M, and the UV luminosity at 1600A, 
Li 600 - To this aim we use the abundance matching tech¬ 


nique (Yuc et al. 2015 Conroy & Wechslor 2009 Behroozi 
et al.|20l0 Vale fc Qstriker|2004 1, which makes the implicit 

assumption that such relation is monotonic. The UV lumi¬ 
nosity function has been measured recently at z > 8 using 
the HST legacy fields ( Bouwens et al.|2014a|b l 
fitted by a Schechter parametrization ( |Schechter|1976| ): 

dn(Mi 6 oo, z 


and is well 


dMieoo 


= 0A\n(10)(p*x 1+o ‘e x , 


(1) 


where x — io°- 4 ( m *6oo M i6oo) an( j M 1600 is the absolute dust 
attenuated AB magnitude at 1600A. Bouwens et al. (2014b I 


computed the redshift evolution of the three Schechter pa¬ 
rameters (cp *, M* 600 , a) with a linear fit of the observational 
data from z = 4 — 8 dropout galaxies, finding significant 
evolution in <)>* and a: 

M* 60 o = (-20.95 ± 0.10) + (0.01 ± 0.06)(z - 6) (2) 

<p* = (0.47lo.iJ)10 (_0 - 27±0 05)(2:_6) 10“ 3 Mpc“ 3 (3) 

a = (-1.87 ± 0.05) + (-0.10 ± 0.03)(z - 6). (4) 

We extrapolated these fits to higher redshifts. With this in¬ 
formation and the mass function dn/dM (Press & Schechter 


19 74] [Shet h & Torrnen |l999 Sheth et al.||2001) we can de¬ 

rive a M 1600 = Mieoo(M, z) by forcing the equality between 
the number density of the most luminous galaxies and most 
massive halos: 


f 


dn(M ', z) 
dM' 


dM' = 




dn{M 1600 ,z) w 

-aM i6oo- 


m 1600 (m,z) dM i6oo 


( 5 ) 

Since the SFR is correlated with the unobscured UV 
magnitude M( 600 , we have to correct for dust. We will use 
the standard correction using the spectral slope pt ( f\ cc X 13 ) 
fitted by Bouwens et al. (2014b I 


P(Mieoo, z) = ,d_ig. 5 ( 0 ) + 


d/3{z) 


Planck Collaboration et al. (2015)). 


(Mieoo + 19.5) (6) 

( 7 ) 

= —0.18 — 0.03(z — 6). (8) 


dMieoo 

0 - i 9 .5(2) = —1-97 — 0.06(z — 6) 
dfijz) 
dMieoo 
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Then we have (Meurer et al. 1999) 


A/l600 = Ml600 — ^1600 
Tieoo = 4.43 + 1.990. 


( 9 ) 

( 10 ) 


The top panel of Fig.[l]shows the resulting Li60o(M, z), 
the intrinsic spectral luminosity density at 1600A, at a = 4, 
6 and 8; the bottom panel shows the reddening E(B-V) (see 
eq. (24)), along with the average weighted with the Lya lu¬ 
minosity (eq. ©)■ At higher redshift galaxies are more ef¬ 
ficient in forming stars and less dusty. 

Using the stellar synthesis code starburst9lf (Lei- 


therer et al.|1999 Vazquez & Lcithcrcr 2005;|Leitherer et al 


20141 we can compute the intrinsic emission of ionizing pho¬ 


tons given the dust corrected UV luminosity: 

r T 1 Mq /yr , 

-Ul600 — -^1600 ^ Y 


finally yielding 


N 912 = N^ 2 @/yl X Vh 


_ m 1 M o/yr 


N 912 = N, 


r lAfg/yr ’ 
^1600 


( 11 ) 

( 12 ) 

(13) 


where V 912 is the rate of emission of ionizing photons. This 
treatment introduces two new parameters, such as the age 
of the stellar population, £ age and the metallicity Z. 

Fig . [2]shows how the ratio L 912 /L 1600 depends on these 
parameters for a galaxy with constant SFR. If we consider 
t agB ~ 10 8 yr, i.e. about 10% of the Hubble time ( th ) during 
the EoR, the emission of ionizing photons depends weakly 
on Z and a factor 2 in t age introduces a order 10% variation. 
Therefore in the following we will simply consider t age (z) = 
O.lTif(z), Z = 0.0004 and a Salpeter Initial Mass Function 
(IMF, Salpeter (19551) between 0.1 and IOOMq, because 
further complications (such as the introduction of a mass¬ 
age relation or different star formation histories) would be 
ill-constrained and would not make our conclusions more 
solid. At z > 10 we imposed t a ge(.z) = 5 x 10 7 yr, in order 
to avoid strong variations in A 912 due to extremely young 
stellar ages. 

We need to set a minimum mass M m i n for star-forming 
halos; this choice is tricky, as feedback effects on the sup¬ 
pression of star formation in small halos are not fully un¬ 
derstood. Although the value of AL m j n is irrelevant at z < 8, 
it becomes fairly important at high redshift (the difference 
between models with M m j n = 10 8 Mq and M m j n = 10 9 Mq 
is about 50% at z = 12), where the bulk of star formation 
takes place in small halos (Yue et al. |2014 ). 

We set 

1/2 /-I , \ - 3/2 


M n 


= 10 M 0 ft 




fim(z) 187T 5 


1 + Z 
10 


(14) 


as our fiducial choiccj^] such minimum mass corresponds to 
halos with virial temperature of 10 4 A, the minimum tem¬ 
perature necessary for hydrogen atomic cooling to become 
efficient. 


1 http://www.stsci.edu/science/starburst99/docs/default.htm 

2 Note that this choice leads to a star formation rate density 
(SFRD), p r , higher than typically found by observations (e.g. 
Bouwens et al. (2014b|) because we are extrapolating the lumi¬ 
nosity function to magnitudes well below the observational limits. 



M[M 0 ] 


Figure 1. Top panel: relation between UV luminosity at 
1600A and halo mass, L 16 OO (Af, z), obtained from abundance 
matching and the HST LBGs from Bo uwens et al.| < |2014bj ). The 
dotted, dot-dashed and dashed lines correspond to z = 4, 6, and 
8. Bottom panel: dust reddening, E(B-V), from eq. § and 
( |24| ). We also show the mean (E(B — V)) /. Ly/i (horizontal lines) 
weighted with Lycr luminosity (eq. (|25|) ) ■ 



Figure 2. Lg 1 2 /L] 600 as function of the stellar age t age and 
metallicity Z. A constant SFR and a Salpeter IMF between 0.1 
and 100 Mg has been assumed. 


2.1.1 Population III stars 


Metal-free (usually defined as PopIII) stars, depending on 
their IMF, can produce more UV photons per baryon with 
respect to standard PopII stars. Since the Lya intensity is 
proportional to the integrated emission of these photons, 
this could ease the detection from high redshift. 

In order to estimate the PopIII emission, we used stel¬ 


lar models from Schaerer (2002). Under the constant SFR 


hypothesis we can compute the number of photons emitted 
per baryon (Pullen et al.||2014): 


N =^ 


J dmf(m)Q(r. 


)r(m), 


(15) 


where: mb is the mean mass of the baryons; f(m) is the IMF; 
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^ 912/^912 

Salp. 

Chabll 

Chablll 

Salp. 

1.35 

3.11 

1.93 

Chabll 

0.75 

1.72 

1.07 

Chablll 

1.20 

2.76 

1.71 

^co'nt/^cont 

Salp. 

Chabll 

Chablll 

Salp. 

1.14 

2.80 

1.53 

Chabll 

0.57 

1.39 

0.76 

Chablll 

1.06 

2.61 

1.43 


Table 1. The tables show the ratio of the emission rate of ionizing 
(top) and E=[10.2, 13.6] UV (bottom) photons between PopII and 
PopIII star. Each row(column) shows a different IMF for PopIII 
(PopII) stars. Salp. stands for a Salpeter IMF ( |SalpeteT]|1955[ ) 
between 1 and 100 Mq; Chabll (Chablll) is the PopII (PopIII) 
IMF from ||Chabrie~ 2003[). 


m* is the first moment of /(m); Q(m) is the emission rate of 
photons of a star of mass m while r(m) is its the lifetime. By 
comparing the emission efficiency of PopII and PopIII stars 
it is possible to simply understand how the Lya intensity 
scales with the different stellar populations. 

Schaerer ( 2002j ) computed Q(m) and r(m) for PopII 


(Z = Zq/ 50), for PopIII stars and for ionizing or Lyman- 
Werner (LW, E £ [11.2,13.6]eV) photons. In order to 
roughly estimate the emission rate of E = [10.2,13.6]eV 
photons, we added a factor 1.4 to Qlw (but this extrapola¬ 
tion does not affect the results of this Section). We consid¬ 
ered also different IMFs, using the PopII and PopIII IMF 
from Chabrier (20031. 

Table [ij shows the ratio of the number of ionizing and 
continuum photons emitted per baryon by PopII and PopIII 
stars. For all the combinations considered, the number of 
photons emitted per baryon is of the same order of mag¬ 
nitude: therefore it is unlikely that PopIII stars could be 
exploited for the detection of the Lya signal at high red- 
shifts. For this reason in the following we will neglect the 
possible presence of PopIII stars in our sources. 


2.2 Escape fraction 


Ionizing radiation is strongly absorbed by H i and dust in 
the ISM. As a result, only a fraction f esc (M, z) of the pho¬ 
tons with energy > 1 Ryd produced in a halo of mass M at 
redshift s manage to escape into the IGM and contribute to 
reionization. 

Theoretical estimates of / eac are difficult, as such quan¬ 
tity depends on complex details of the ISM structure in 
galaxies. At low redshift (z < 4) / eac is experimentally con- 


strained to be < 0.1 (Steidel et al.|2001 

Shapley et al. 2006 

Inoue et al. 2006; Siana et al. 

20071 llwata et al. 2009; Bout-| 

sia et al. 2011:; Vanzella et al. 

2012a b). 

As we will see later, 


in this redshift range the exact value is not very important 
because the results only depend on (1 — / esc ) s=s 1 , but dur¬ 
ing EoR. different values of / eac lead to considerably different 
reionization histories. 

To further complicate the above situation, results from 


numerical simulations are often inconsistent 1 

Dove et al. 

2000; [Wood & Loeb 2000; Gnedin et al. 2008 

Razoumov 

& Sommer-Larsen 

2006 

Wise & Cen 2009 

1, even if they 


& Loeb (2013) suggest a bimodal function, with / ei 


for the smallest galaxies and / eac ~ 0 for the others. In this 
work we will use this approach, tweaking the three free pa¬ 
rameters (the two escape fractions for small and big galaxies 
and the threshold mass) to have a reasonable reionization 
history: 


f 0.05 for M > 2 x 10 9 M 0 
\ 0.35 for M < 2 x 10 9 M 0 


(16) 


Therefore the effective ionizing rate measured in the IGM 
from a given halo is / eac (A/, z)N 9 i 2 (M, z); the density rate 
of ionizing photons is then obtained from 

h ion = f dMf esc (M,z)N 912 (M,z)-^-, (17) 

The mass-weighted (/ esc (z)).M approaches 0.35 at very high 
redshift because the mean mass of the halos decreases be¬ 
low the threshold mass. As we are going to see next, this 
behavior accelerates reionization and suppresses Lya emis¬ 
sion from low-mass galaxies. 


3 REIONIZATION 

A detailed treatment of reionization is not essential in our 
analysis, because, as we will see in §4.2| recombinations from 
the ionized IGM are not the dominant source of Lya pho¬ 
tons. It could appear that the computation of the neutral 
fraction, xhi, in the IGM is essential, because Lya photons 
interact strongly with H i . However, although numerical 
simulations indicate that H i can reduce significantly the ob- 


served luminosity of individual galaxies ( 

Dayal et al. 

2011 

Mesinger & Furlanetto|2008| |McQuinn et 

al.|2007j Dijkstra| 

2014 

1, this effect is not relevant for Lya intensity mapping. 


This is because Lya photons are not destroyed by H i , but 
are simply scattered; the net effect is a damping of fluctu¬ 
ations on scales below the typical photon mean free path, 
A « 1 Mpc. 

The diffuse IGM emission depends on the opacity of res¬ 
onant Lyman photons which is derived from observations at 
a < 5.5 and extrapolated at higher redshifts. We use a sim¬ 
plified approach to compute the ionization state of the IGM 
that assumes a bimodal (two-phase) distribution in which 
the gas is either fully ionized or neutral; in other words, we 
ignore for the moment the residual H i in the ionized regions. 
We will further discuss this issue in Sec. 14.31 

The purpose of this Section is to compute the filling 
factor Q, e.g. the volume fraction filled by ionized IGM. 
The governing evolutionary equation can be obtained by 
imposing a detailed balance condition among the number 
density of ionizing photons, ionizations and recombinations 
( [Loeb fe Furlanetto|[20l3| : 

nH ~^~ = ™ ion ~ n Hp + z) 3 C(z)a B (T)Q, (18) 


where C(z) is the fiducial clumping factor computed by|Paw- 


lik et al. (2009) at z > 6 (it was used the fit with Cioo 


and z r = 9.0 from appendix A, with C(z < 6) = C(6)) 
and as(T) is the Case B recombination coefficient ( |Meiksin| 
2009:; Loeb & Furlanetto 2013 ). Fig. [3] shows the numer¬ 


ical solution of eq. (181. We used the reionization history 
model to constrain the escape fraction; we obtained Q = 1 at 
z rts 6.0, Q = 0.9 at z « 6.4 and Q = 0.5 at z = 8.6. The cor- 
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Figure 3. Redshift evolution of the neutral hydrogen fraction 
2?hi = n Hi/ n H (dot-dashed line) and residual x^u. in the ion¬ 
ized regions (solid line). During reionization a?Hi ~ 1 — Q; after 
its completion, x^i is computed from the ionization equilibrium 
condition matched to quasar absorption line data (see eq. 1 35 [ ) . 


Figure 4. Average Lya escape fraction, (fe SC )Lya from |Atek| 
|et al.| < |2009[ ) (orange dash-dotted, eq. ( |23| >) weighted with 
the Lya luminosity (eq. |25 >); for comparison, the average 
Calzetti dust attenuation law (Calzetti et al. ]2000) (eq. (gl}) at 
912 A(black dashed) and 1216 A (black solid). 


responding electron scattering optical depth is r es = 0.067, 
compatible with the most recent Planck result ( |Planck Col-| 
laboration et al.|2015 1 of 0.066 ± 0.016. 


4 MEAN Lya INTENSITY 


In this Section we will compute the mean intensity of the 
Lya signal. The intensity, /, observed at frequency Vo and 
redshift a = zo is simply the integral of the emissivity, e, 
along the cosmological path £ (Loeb & Furlanctto 20131: 


I{v 0 , z 0 ) = L(i + Zo ) 3 J* dz'^e [i/ 0 ( 1 + z),z] 


(19) 


where d£/dz = c[(l + z)H{z)\~ 1 . In this context the exact 
shape of Lya line profile is not relevant as in practice the 
typical line width is negligible compared to the spectroscopic 
resolution of the instruments. For this reason we assume a 
Dirac delta-function, Sn- Hence 


z) = p Lva {z)5 D {v - u a ), (20) 


p Lya is the Lya luminosity density. In this case the integral 
reduces to 




c p Lyoi (z a ) 
4tt v a H(z a ) ’ 


( 21 ) 


where z a = Va/vo — 1 is the redshift of emission of the 
Lya photons. We will consider three main processes lead¬ 
ing to Lya emission: recombinations in the ISM of ha¬ 
los/galaxies (Phaio)i recombinations in the IGM (Piqm)> an d 
excitations due to UV photons (p^on “). 


4.1 Recombinations in the ISM 


from the 2 P 2 state a Lya photon is emitted; (iii) decay from 
higher levels produces the emission of a Lyn photon, or (iv) 
a cascade to another excited state and then one of the first 
three cases. Lyman photons interact strongly with hydrogen 
and we can assume the Lyn photons to be locally absorbed 
and reemitted as long as they end up into a Lya photon or 
in a two-plioton emission. If we include all the four emission 
mechanisms, a fraction f a of the Case B recombination emits 
a Lya photon ( |Canta lupo et al. 2008]) 

f a = 0.686 - 0.106 log T 4 - 0.009TT 0 - 44 , (22) 


where T 4 = T/(10 4 K). As the temperature dependence is 
weak, we will assume in the following T 4 = 2 and hence 
f Lya = 0.65. 


As we have seen, a fraction f esc (M,z) of the ioniz¬ 
ing photons escape into the IGM, while the remaining 
part are absorbed by H 1 photoionization. The recombina¬ 
tion timescale for typical ISM densities (uh ~ 1 cm -3 ) is 
t rec S3 l/niias « 2 x 10 5 yr and is negligible compared 
to the Hubble time. Thus we can assume photoionization 
equilibrium and that any H i photoionization generates f a 
Lya photons. 

However, the effects of dust absorption can be very sig¬ 
nificant especially at low redshift. Depending on the dumpi¬ 
ness of the ISM, the resonant nature of the Lya line can lead 
to both a suppression ( Hayes et al.|2011 1 or an amplification 
( Neufeld|1991 Verhamme et al.|2006 1 of Lya radiation with 
respect to the UV continuum. To overcome these difficulties 
we resort to an empirical approach. Atek et al. (2009) found 
a relation between the Lya escape fraction (/“ c ) and the 
reddening E(B — V) in low redshift LAE (z tts 0.3). 


= w -o.ik Lya E(B-v), 


(23) 


When an electron recombines with a proton, it can directly 
populate the ground state (Case A) or an excited state (Case 
B); in the second case there is a chance of Lya emission. 
Basically, an excited hyd roge n atom can decay to the ground 
state in four ways (e.g. Dijkstra 2014)): (i) from the 2 S 2 
state, it can decay only through two-photon emission; (ii) 


where kL ya = 12.7. We computed the E(B — V) of each halo 
from the observational /3 (eq. ©) and using the relation 
(Calzetti et al. J2000) 

E(B - V)(M, z) = 0.44Ai 6OO (M, z)/ 4.39. (24) 

As a caveat we point out that the we are implicitly making 
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two important assumptions: (i) f^ sc depends only on the dust 
content, (ii) we can extrapolate results obtained for low- 
redshift systems to high-z. In spite of these uncertainties, 
the present approach has the advantage that it yields the 
extinction of Lya photons by dust as a function of mass and 
consistently with the UV dust absorption. Fig. 4 shows the 
average Lya escape fraction, (fesc) Lya- 

Under these hypothesis, the halo luminosity is given by 

LLya (A/, z) — hpULyafa^ 

/“ c (M,z)( 1 - /esc(M,z))JV 912 (M,z); (25) 

the luminosity density is 

&= [ dM^L Lya (M,z), (26) 

J ^min 

And, finally, the intensity is 

few - OM^d-MN n , (27) 

Clearly, considerable uncertainties remains on the parame¬ 
ters in this formula: we have already discussed those con¬ 
cerning fe SC and / esc ; we also recall here that Ngi 2 is ex¬ 
trapolated to masses not probed by HST surveys. These 
uncertainties are amplified in the estimate of the PS, which 
scales as iLya 2 - Therefore, from this estimate we cannot 
safely predict the Lya PS that should be observed, but we 
can conclude that a survey can constrain the mass-averaged 
values of /“ c , / es c and IV 912 . Of course, these parameters are 
strongly degenerate. 


4-1.1 Lya luminosity function 

As a by-product of the method and a sanity check, we can 
predict the Lya luminosity function (LF) of our galaxies. 
This is shown in Fig. [ 5 ] for z « 5.7 (solid line), where it is also 
compared to the experimental data from Ouchi et al. (20081; 


although the agreement seems adequate, we recall that we 
have not included the Lya line damping (Dayal et al.|20lT 


Dijkstra|2014 1 of the blue part of the line due to the inter¬ 
vening IGM. To provide a rough estimate of this effect we 
also show the LF in which a constant transmissivity factor 
T a = 0.5 (dotted line). Additionally, we note that if LAE 
are characterized by a duty cycle, this could also modify the 
shape of the LF at fixed p^aio ■ For example the dashed line 
in Fig. [ 5 ] corresponds to a LF with: (1) T a = 0.5, (2) a 30% 
duty cycle (i.e. only 30% of the galaxies host a LAE), and 
(3) an increased (by a factor 3.3) individual LAE luminosity. 
Using more sophisticated (but still arbitrary) prescriptions 
for the duty cycle it should be possible to obtain a LF in per¬ 
fect agreement with the observational data. For the present 
goals, we consider the accuracy of the present results as sat¬ 
isfactory. 


4.2 Recombinations in the IGM 

In our reionization model, we approximate the IGM in 
a two-phase medium defined by its ionization stata (neu¬ 
tral/ionized). Recombinations occur only in the ionized re¬ 
gions, and produce a luminosity density: 

Pigm = h p v a faQ{z)n 2 H {Q){l + z) 3 C{z)a B (T); (28) 



L[e rgs J ] 

Figure 5. Lya luminosity function at z = 5.7. We show our 
model (i) without IGM absorption or duty cycle (solid line); (ii) 
with a 50% IGM absorption (dotted line); (iii) with also a 30% 
duty cycle (dashed line). The experimental data are from [Ouchi| 
|et a l. ( 2008[ l. The faintest galaxies in our model have luminosity 
ps 9 X 10 3 ^erg/s and therefore are not shown in this plot. 


note the different dependence from the ionized fraction com¬ 
pared to the emissio n from a uniformly ionized gas which is 
oc Q 2 (z ) ( Pullen et al.|2014 i. 

It follows that the intensity is 


Pgm( z ) — 


chp 
47 vH(z 


-faQ(z)n H (0)n H (z)C(z)a B (T). (29) 


4.3 UV excitations in the IGM 

Continuum UV photons with energy 10.2 — 13.6 eV (i.e. be¬ 
tween the Lya line and the Lyman-limit) after escaping from 
galaxies redshift into Lyman-series frequencies and can be 
absorbed by H 1 . The resulting excited H 1 atoms eventually 
decay to the ground state, producing Lya photons. The com¬ 
putation of the Lya emissivity in this case is slightly more 
involved than for recombinations, as it requires a model for 
the production and absorption of 10.2 — 13.6 eV photons 
( Pullen et al.|2014 |. 

The absorption probability of Lyn photons at redshift z, 
P a bs(n, z), can be expressed by the canonical Gunn-Peterson 
( Gunn &; Peterson|1965 I optical depth: 

r G p(n,z) = ^jrz\nm(z)\ 3 Lyn , 


2 H(z) 


where 


7 ™ = 50 


1 — n 
0.75 


fn 


0.4162 


MHz 


(30) 


( 31 ) 


is the resonance half width at half-maximum, f n is the os¬ 


cillator strength tabulated in Verner et al. (1994), and A l 


is the wavelength corresponding to the n-transition. 

A precise computation of top requires the knowledge 
of riHi even in ionized regions. Our reionization model is 
too simple to compute it; but quasar absorption line obser¬ 
vations provide the “effective” optical depth (i.e. the aver¬ 
age of tgp over small scale inhomogeneities of overdensity 
A = p/p) for Lya photons up to z = 5.5 (Fan et al. 2006) 

4 . 3 ± 0.3 


Tefffc = <TGP (2, z)> A = (0.85 ± 0.06) 


IT 2 


(32) 
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z 

A 

P 

Co 

2 

0.170 

2.23 

0.623 

3 

0.244 

2.35 

0.583 

4 

0.322 

2.48 

0.578 

6 

0.380 

2.50 

0.868 

8 

0.463 

2.50 

0.964 

10 

0.560 

2.50 

0.993 

12 

0.664 

2.50 

1.001 

15 

0.823 

2.50 

1.002 


Table 2. Parameters for the IGM overdensity probability distri¬ 
bution function, P( A, z), used in eq. ( ].‘{4| : data are from |Miralda-| 
|Escude et al.| (|2000[) 


For z > 5.5 even the ionized IGM is optically thick to the 
Lya line, which is the most relevant excitation (in the spec¬ 
tra from starburst99 the bulk of redshifted UV photons 
has a frequency between the Ly/3 and the Lya one). Our 
strategy is then to use a model for small scale inhomo¬ 
geneities from which we compute Teff, a as a function of nm 
in ionized regions. Then nm is tuned to reproduce 
The extrapolation to higher redshifts is arbitrary but not 
relevant for the results of this work, because the absorption 
is already saturated at a = 5.5 in the Lya line. 

To understand how tgp depends on IGM (over) density 


A we impose photoionization equilibrium (Choudhury & 


Ferrara 2005), an excellent assumption under IGM condi¬ 


tions: 


uhiT = nmin e as(T)C'(z), (33) 

with nm = (1 — x e )riH, n H n ~ n e = x e n,H , a B {T ) oc T -0 '. 

As we are concerned with ionized regions in the post¬ 
reionization epoch, x e ~ 1 ; in addition, we use an equation 
of state T oc A 7-1 , where 1 < 7 < 2.4 is the adiabatic 
index. The results depend weakly on 7 and we assume an 
equation of state with 7 = 1.4. With these assumptions, 
tgp oc A 2 -' -0 ' 77 , i.e. tgp oc A 1 ' 2 for 7 = 1.4. 

As a model for the small scale inhomogeneities of the 
IGM, we adopt the one proposed by |Miralda-Escude et al (| 
(|2000). Such model, tested against numerical simulations, 
gives the Probability Distribution Function for the overden¬ 
sity A as 


P( A, z) = A(a)exp 


[A 2 / 3 -C 0 (z)] 2 \ p 
2[2A 0 («)/3 ] 2 ) 


(34) 


where the constants A, /?, Go are reported in Table [2] as a 
function of redshift and 5o(^) = 7.61/(1 + z). 

The effective optical depth can now be computed as 


g—effKD = ( e -^ A = / p(A,z)e- Al ' 72 TGpM dA. 

Jo 

(35) 

We have set Ah = 100 as a reasonable overdensity of the 
external parts of virialized halos; the results are however 
insensitive to the precise choice of this parameter. Using 
the least square method and eq. ( |32[ ) between z = 4 and 
5.5 we find xm = 1.1 x 10 _5 [(1 + z)/ 5 ] 4,53 (Fig. [|]) and 
extrapolate this relation to z > 5.5. Given that the bulk of 
the continuum photons is emitted between the Lya and the 
Ly/3 line and that the absorption probability of Lya photons, 
P 7 b *“, saturates at 2 = 5.5, the final results do depend only 
weakly on this extrapolation. 


K! 


Figure 6. Transmissivity of different Lyman-series lines through 
the IGM as a function of redshift. During the EoR the Lya line 
saturates and is completely absorbed. 



If we consider the neutral IGM ( x e C 1), we have in¬ 
stead tgp oc riHi oc A, or 

e -Teff,neut(",*) = f ^ p( A , z ) e - Atgp (n ’ z) dA; (36) 

Jo 

this results in P abs = 1 for the lowest Lyman levels. 

Finally, the transmissivity of a Lyn photon through the 
IGM at redshift z, T(n,z) = 1 — P abs (n, z), is: 

T(n, z) = Qe~ T ‘ i!M + (1 - Q)e~ T ^ (37) 

Fig .[ 6 ] shows the above result in graphical form for different 
Lyman-series lines. It is evident that during EoR the IGM is 
optically thick at least to Lya photons, therefore all contin¬ 
uum photons are scattered and re-emitted as Lya photons. 
On the other hand, the increasing transparency of the IGM 
at z 6 — 7, where T(n, z) approaches unity, has important 
consequences that were not properly considered in earlier 
works (Pullen et al. 2014). 

A fraction /(n) of the excited H i atoms produce Lya 
photons (see Table 1 from Hirata (2006)), therefore the Lya 
intensity at redshift a is 

OO U 

I{u Lya ,z) = V lpVLva P abs (n, z)f(n)I(v Lyn ,z); (38) 
7^2 bpVLyn 


Then, from the radiative transfer equation ( Loeb fe Furlan-| 
|etto|2013| ) the intensity observed at z — 0 and v = vt , ya /(1 + 
z) is 


rLya 


{z) = I(v Lya ,z)/(l + zf. 


(39) 


Hence, the final step is to compute I{uL V n, z); this is 
easily done once the UV emissivity, e(p), is assigned. We 
will parametrize the specific emissivity as: 


e(i/, z) = hpvh v (v, z) = 


-/ 

J M 


dM J^M M ,z)f^Zt(M,z), (40) 


where L„(M,z), the spectral luminosity density, is com¬ 
puted with starburst99 coherently with Sec. |2.1| and 


/dust is the fraction of UV photons that escape unimpeded 
by dust. The latter can be accounted by using, e.g. the 
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Calzetti extinction law ( Calzetti et al.pOOO Li et al.||2008 l 
in parametric form: 




(41) 


where E(B — V)(z) is given in eq. \2A\ . The redshift evo¬ 
lution of (fdust{ z ))Lya is plotted in Fig. [ 4 ] for non-resonant 
1216 A photons and for 912 A ones. In this paper there are 
not significant differences between the dust attenuation of 
resonant Lya photons and the neighbouring ones. 


For the calculation of the intensity, in this case we have 
also to consider absorption from H 1 ; indeed a photon, while 
redshifting into a Lyn frequency, can pass through other 
Lym frequencies, with m > n. Therefore we will include a 
product of the transmission fractions from eq. (|37|): 


I(yLym z ') — tlpVLyn 


P + z ? 


47T 


f 




dz' x 


ch„(v', z') 
H(z') 


n 

ri/=n+1 


') 


T(n' ,z n i(v' ,z')), (42) 


where (i) v' = VL yn {l + «')/( 1 + z) is the frequency 
of a photon at z' that redshifts to vp yn at z; (ii) 1 + 
Zmax(n,z) = {v L L/vLyn)(l + z) is the redshift at which 
a Lyman limit photon redshifts into a Lyn photon at z; 
(hi) n' max = [1 - (u Lyn /u LL )(l + + z)}~ 1/2 is the 

maximum m such that v’ > i/Lym, and (iv) 1 4 - z n / = 
(I'Lyn 1 / V ')P + z ') is the redshift at which a photon with fre¬ 
quency v' from redshift z' redshifts into the Lyn' frequency. 
Thus from eq. ([39|l: 


^out i.z) = i F Pabs(n,z)f{n)x 


/ 


dz 


,n v (y', z') 
H(z') 


n 

n , =n -\-1 


(43) 


As in Sec. 4.1 in eq. {431 there are several poorly constrained 


parameters, both in n„ and in the IGM modelling. There¬ 
fore, if UV excitation dominates the Lya PS, its measure¬ 
ment can constrain the intrinsic emission of UV photons, 
dust obscuration and the ionization state of the IGM. As 
for the halo emission, these three physical processes are de¬ 
generate and cannot be disentangled. 


As a final remark, we note a subtle point concerning 
photons scattered directly by the Lya line: if the IGM were 
perfectly transparent, we would observe exactly the same 
intensity at the redshifted Lya frequency. The effect of the 
IGM is to substitute part of the resonant Lya photons from 
the line of sight (LOS) with the local average over the solid 
angle: 

A q (IIlOs) —> Pabs(A c ,)n + (1 — Pnbs)Iv a (^LOs)- (44) 

In a perfectly homogeneous universe these two terms cancel 
each other out; however, for spatial fluctuations the can¬ 
cellation does not occur and also the continuum emission 
along the LOS is suppressed by the foreground removal. 
Fig-0 shows the difference in Lya intensity when the IGM 
Lya scattering is either included or excluded. 


4.4 Additional processes 


Atoms can emit Lya photons also when they are collision- 
ally excited and emit decaying to the ground state. This 
mechanism needs both thermal kinetic energies comparable 
to the Lya energy and the presence of H i . Such condi¬ 
tions can be reasonably found only in the transition zones 
between ionized and neutral regions. Thus, the emissivity 
depends strongly on the morphology and temperature of 
the ionized regions, which is hard to model analytically. 
Moreover, as the excitation temperature of the Lya line, 
T a = 1.2 x 10 5 K is much larger than those typical of photo- 
heated gas (~ 1 — 5 x 10 4 K), the abundance of energetic elec¬ 
trons is exponentially suppressed. |Cantalupo et ah] ( |2008| ) 
computed the collisional excitation coefficient for Lya emis¬ 
sion (along with the equivalent one for recombinations 

^Lya — f Lya^B ) 


eff . eff 

riLya = n e nmia Lya + n e nmqLya- 


(45) 


As expected, depends strongly on temperature: at 

T = 2 x 10 4 K it is three orders of magnitude larger than 
o^Lya- Therefore at this temperature collisional excitations 
are relevant only if xhi A, 1CU 3 ; however, a factor 2 change 
in temperature results in a variatio n of three orders of mag- 

12013) assumed the collisional 


nitude in q e A 


Silva et al. 


excitation intensity to be 10% of that from recombinations; 


Pullen et al. (2014) found it negligible and did not consider 


it in their model. To avoid the introduction of poorly con¬ 
strained parameters related to the morphology and the tem¬ 
perature of the gas, we make the most conservative choice 
and neglect this contribution. 

An additional ionization source could be represented 
by exotic annihilating/decaying dark matter particles, that 
partially ionize the IGM almost uniformly. In this case the 
additional Lya intensity from the IGM is 


risa/,', 1 Q( z ) _ / . \2 

uni K > O(z) 


L X( Z Y 


2 x 10 4 A' 

T...: 


0.7 


■tlGM,: 


c(*). 


(46) 

The temperature dependence (T un i) comes from the 
recombination coefficient as and introduces a potentially 
large factor for a luke-warm gas heated by X-ray photons; 
x(z) is the ionized fraction outside the reionized bubbles. 
Given the uncertainties in the nature of the dark matter 
particle and annihilation cross-sections, we will not consider 
this source in detail; however, we point out that it might be 
significant before the completion of reionization (Fig. 0 ). 


4.5 Total intensity 

Fig- 0 show the results of the previous estimates. The dom¬ 
inant source is represented by the scattering of continuum 
photons. This can be understood from the fact that stars 
emit more UV photons with 10.2eV < E < 13.6eV than 
ionizing photons (see Fig.0. The efficiency of these UV pho¬ 
tons in producing Lya photons is comparable to the ionizing 
ones and the effect of dust absorption is of the same order 
of magnitude. However, as we will see in the next Section, 
fluctuations induced by LTV photon scattering have smaller 
amplitudes on small scales, because while UV photons origi¬ 
nate from galaxies, Lya photons are eventually produced in 
more distant regions of the IGM. 
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Figure 7. Lya mean intensities. Each curve represents the con¬ 
tribution from different physical processes: (i) recombinations in 
the ISM (solid), (ii) recombinations in the IGM (dot-dashed), (iii) 
UV excitations in the IGM (thick dashed), (iv) same as (iii) but 
without photons scattered in the Lya line (thin dashed), (v) emis¬ 
sion from an hypothetical uniformly ionized cold IGM (xhi = 0.3, 
Tjgm = 100K, dotted). 



Figure 8. Ratio between the emission rates of E = 10.2 — 13.6 
eV and ionizing photons. A constant SFR and a Salpeter IMF 
between 0.1 and 100 Mq has been assumed. 


Recently Croft et al. (2015) made the first intensity 
mapping observation by cross-correlating the spectra from 
the SDSS with the QSO distribution at 2 = 2.55. They 
found a Lya comoving luminosity density of Pgnss ~ 3.1 x 
10 41 (3/6 a )erg s _1 Mpc _1 (b a = 3 is their fiducial mean 
Lya bias), with statistical uncertainties of ~ 25%. Our 
model predicts Piqm = 8.73(6.51) x 10 4O erg s _1 Mpc _1 and 
Pism = 6.62(3.21) x 10 40 erg s _1 Mpc _1 at z = 4(7). As we 
will see in Sec. |5.4[ in this context only the ISM emis¬ 
sion is important, since the IGM emission is smooth on 
the scales considered by |Croft et al.| ( [2015] ) (< 100 Mpc). 
If we linearly extrapolate our results to z = 2.55, we find 
Pism(z = 2.55) « 0.91 x 10 41 erg s _1 Mpc _1 , a value of the 
same order of magnitude of the Croft et al.| (j2015| results, 
albeit a factor of 3 smaller. 


5 INTENSITY FLUCTUATIONS 


The observation of the mean intensity of Lya emission is 
severely limited by foregrounds: the intensity of the NIR 
foreground around lpm is of about 10 3_4 nWm _2 sr _1 
IKashli nsky|2 0051, while our predicted Lya signal is 5-6 or¬ 
ders of magnitude fainter. On the other side the situation is 
different if we consider the fluctuations, because the dom¬ 
inant foreground components are smooth in frequency and 


can be removed (Wang et al. 2006 Chapman et al. 2015 
Wolz et al.|2015 i. 


In this Section we will develop the theory necessary 
to predict the power spectrum of the fluctuations with our 
model. As for the intensity we will examine the three main 
Lya photon emission processes (ISM emission, IGM recom¬ 
binations and IGM Lyman absorption) separately. 


5.1 Fluctuations from ISM recombinations 


The distribution of galaxies in space fluctuates due to: (i) 
clustering of large scale structure, and (ii) shot noise due to 
the discrete nature of the sources. 

The clustering term can be written in terms of the halo 
mass function. In the linear regime, 5(dn/dM) ss b(dn/dM)5 
(Sheth & Tormen|l999 Sheth et al.||200jj . Thus, 

57 haio(«> x ) = T hZ 0)( 6 (z)M(x), (47) 

where (b(z)) a is the mean bias weighted with the Lya lumi¬ 
nosity and the mass function (Fig. [9|: 


{b(z)) a = 


J+2 Ii ^ ya (M,z)^b(M,z) 

• Lya 
Phalo 


(48) 


For simplicity we have neglected redshift space distortions, 


because the interplay between peculiar velocities (Kaiser 
19871 and radiative transfer effects (Zheng et al.|2011 Croft 
et a.1.] [2015]) is unclear. This is a conservative choice be¬ 
cause this effect increases the linear bias (Gong et al.||2014 


Visbal et al. 20151. Therefore the Lya power spectrum 


from galaxies is proportional to the linear power spectrum, 


P gg = (^)a cr fc ( we used the transfer function from Eisenstein 


& Hu 19981). It follows that 


= &z) 2 m)la 2 k (z). 


(49) 


We have also neglected fluctuations associated with the dis¬ 
tribution of satellite sources within a single halo (ofter re¬ 
ferred to as the two-halo term) because on those scales 
(k > hMpc -1 ) the power spectrum is completely dominated 
by shot noise. 

The second source of fluctuations, shot noise, is due 
to the Poisson fluctuations of the number of galaxies in a 
volume V around the mean value, 


_. dn , _ , 
N m = V -jrzdM 
dM 


with M £ [M, M + dM], 


(50) 


and with standard deviation ( 5N fj) = = Nm- To com¬ 

pute the power spectrum we can imagine to divide the mass 
range into a large number of bins of width dM centered at 
Mi (i 3> 1 is an integer index). As different bins are uncor¬ 
related, and recalling that 

SNMi 
V ’ 


5 S n(*,x, dM) = 


(51) 
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we can directly write the expression for the luminosity den¬ 
sity fluctuations: 

(52) 

More often, though, we will be interested in the correlation 
function. Starting from 

( 5 sn(*,x, dM)5sa{j,-x!,dM)) = - x')dM, 

(53) 

we then find 


( 5 Phalo(«. X ) 5 Phalo(^> X ')} = 
= PhaTO 


_ L Lya (MQ 2 Jj,(x - P)jFidM 

•Lyot 




(54) 


Finally we take the continuum limit dM —> 0 
to obtain the shot-noise power spectrum. Since 
( s Ph^o( z ’ x ) s Ph^( z ’ x ')) x <&(x-x'), the power spectrum 
is independent of k, for the lack of spatial correlations: 


f 

pLya ( *. _ rLyotf \2Jl 
Ml,SNW — -'halo 


mZ7 d M^L Lya (M,z)\ 

• Lya/ \2 ’ 


•Lwa/ 
PhaloW 


(55) 


Shot noise dominates the signal for k > 1/iMpc - , but 
on these scales the luminosity from the point sources is 
smoothed by the scatterings in the CGM/IGM, suppress¬ 
ing the power spectrum. Thus, while measurements of the 
clustering component on scales k < 1/iMpc^ 1 are unaffected 
by this problem, detections of the Lya signal, as proposed 


by e.g. Pullen et al. (20141, might run into difficulties as the 


bulk of the predicted signal to noise originates from uncer¬ 
tain high-fc fluctuations. It is very difficult to investigate this 
effect analytically, because it depends both on the distribu¬ 
tion of the relic H i in the ionized regions and on the size of 
the ionized bubbles. As we have not modeled these aspects 
in detail, our shot noise results should be used with caution. 


5.2 Fluctuations from IGM recombinations 


IGM emission is powered by hydrogen recombinations; to 
compute its fluctuations we then need to start from the hy¬ 
drogen density field, 8 h■ Although assuming that baryons 
trace dark matter is a common practice, this hypothesis 
becomes increasingly inaccurate on small scales, where gas 
pressure can smooth fluctuations. This effect can be appro¬ 
priately described by a simple model from [McQuinn et al.| 
(2005), prescribing that at 2 /5 8 fluctuations are truncated 


for k > k F (z): 


P H (k, z) = F h {k/k F {z))al 


(56) 


Fb ^ 2 ® + 2(l + 4a; 2 ) 1 /4 ! ( 57 ) 

where k F (z) = 34fl m (2) 1 ^ 2 hMpc _1 and Fb{x) = 1 for z > 8. 
This correction, however, turns out to be irrelevant for our 
work because we will be mainly concentrated on large scale 
structure signales, i.e. fluctuations with k fa 0.1/iMpc -1 , and 
therefore Sh ~ 8 . 

The intensity from IGM recombinations is proportional 
both to the square of the hydrogen density and to the filling 


factor (eq. (291): 


Agm « Q(, z )n H ; 


(58) 


These two terms generate a different types of fluctuations: 
hydrogen density fluctuations are the only ones present after 
the completion of reionization; those induced by the patchy 
reionization morphology dominate at high-redshift, as ion¬ 
ized regions are biased. 

We will first analyze fluctuations in the hydrogen den¬ 
sity. In the linear regime 8n 2 H w 2 8uh and: 

«]££(*, x ) « 2«5(x)/ IGM (59) 

Therefore the power spectrum is simply: 

p igm ,m(*5 z ) = 4 (Agm ) 2 °k( z )- (60) 


The treatment of fluctuations in the ionized fraction is more 
complex and has been discussed by Santos et al. (2003); we 
adapt that method to our purposes. For overdensity S, we 
can write SQ(z, 8 ) fts Q(z)( 1 + 8 q); if further <5 1, then 

5q oc 8 . Therefore, we can define the porosity bias as 


- S Q _ Q( Z ’ S ) , 

W * )= T-W ' 


(61) 


To compute Q(z, 6 ) we use the theory developed in Sec. [3] 
allowing for an explicit dependence on 8 in eq. (181: 


dQ{z,5{z)) _ (1 + (b(z))gi 28 (z)) dQ(z,5 = 0) 
dt riff(0)(l + 8 (z)) dt 

- a B C(z)n H { 0)(1 + <5(«))(1 + z) 3 Q(z, 8 (z)), (62) 


where (b(z ))912 is the mean bias weighted with the number 
of ionizing photons emitted by the halos (Fig. [91 and 8 (z) = 
8 (z = 0)D(z), where D{z) is the growth factor (Liddle et al. 
19961. 


Finally, we have to account of the fact that the relation 
between Sq and 8 is not linear on all scales: if ionized bubbles 
have a typical radius Rb, fluctuations on smaller scales must 

1 2 n / \ 2 

be damped. For this reason we add a factor e~ nb(z> to 
our final power spectrum. Our formalism is not suited to 
compute the typical radius of the ionized regions; thus we 


will use the same crude model of Santos et al. (2003): 


Rb = 


1 — Q 


1/3 


Ra , 


(63) 


where Rq = 100 kpc is a guess for the size of the ionized 
bubbles at the beginning of reionization. The final results 
are not strongly dependent on the detailed form of Rb be¬ 
cause we expect that when the bubbles reach a radius of 
interest in this work (e.g. of the order of tens of Mpc) the 
ionization fraction is almost unity and therefore this source 
of fluctuations is negligible. Moreover this transition is very 
fast and relevant only in a limited redshift range. 

Fig. El shows the porosity bias 6 q as a function of 
redshift along with bQe~ k Rb ^ (k = lMpc^ 1 ). It is evi¬ 
dent that the finite bubble size correction is negligible for 
k < hMpc -1 . In addition, note that is always much 
smaller than the average galaxy bias and that the Lya in¬ 
tensity is fainter. Therefore the power spectrum of IGM 
recombinations is negligible when compared to the other 
Lya sources. 

The final power spectrum arising from fluctuations in 
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Figure 9. Bias of the fluctuations in the filling factor used in eq. 
{64}. We show 6 q (orange dotted) and bQe~ k R b/ 2 (blue solid) 
with k = IMpc -1 . We show the mean galaxy bias ( b) a (red 
dot-dashed) from |Sheth et al.| {2001} weighted with the Lya lu¬ 
minosity (eq. @) and the mean bias weighted with the number 
of ionizing photons emscaping into the IGM ( 6)912 (red solid). 

the ionization field is: 

p iGM, xm ( k ’ z ) = (4m f h Q{zj 2e ~ k2Rb{zr vl(z)-, (64) 

In this case it is not necessary to include the shot-noise term, 
because fluctuations in the IGM density are not associated 
to discrete sources and our model is not well suited for the 
calculation of the shot noise arising from the discreteness of 
the ionized bubbles. Moreover, this source of Lya photons is 
sub-dominant and shot-noise is never relevant on the scales 
of interest for intensity mapping. 


5.3 Fluctuations from IGM UV excitations 


Fluctuations associated with UV photon excitations are 
seeded by various contributions. First we consider the bi¬ 
ased fluctuations in the UV emissivity: 


^ontOu x ) = ^ J dQ. j dz A{z,z')n v {y' ,z')x 

x {{b{z')) u '5(z\x +l(z,z')Cl)], (65) 


where l(z, z') is the comoving distance traveled by a photon 
between redshift z and z'; A(z, z') is a complex function 
that can be recovered from eq. ( |43| ), and {b(z'))^i is the 
mean bias weighted with the luminosity at the appropriate 
UV frequency 


= 


r +00 

fM min 


jL v > ( M , v', z')b(M, z')dM 


hpvh u (y', z') 


( 66 ) 

The integrand in eq. (651 has a non-trivial spatial depen¬ 
dence that tends to suppress small scale fluctuations. In¬ 
deed, the bulk of the photons that redshift into the Lya line 
are emitted between the Lya and the Ly/3 line; but the dis¬ 
tance l(z,z') required for a Ly ft photon to redshift into a 
Lya photon can be quite large (~ 575(460)Mpc at z = 4(7)). 
Therefore if we consider fluctuatons on smaller scales, part 
of the UV photons contribute less to the power spectrum 
because their fluctuations do not correlate. This effect is 



k = 0.01 

0.1 

1 

z = 4 

4.72 x 10 " 6 

9.76 x 10 “ 4 

1.71 x 10 “ 2 

6 

5.50 x lO " 7 

6.50 x 10 “ 5 

1.15 x 10 -3 

7 

2.00 x 10~ 7 

2.09 x 10 " 5 

3.63 x 10 -4 

8 

7.67 x 10-® 

7.52 x 10 " 6 

1.27 x 10 “ 4 

10 

1.35 x 10 “ 8 

1.19 x 10 -6 

1.92 x 10 “ 5 


Table 3. Predicted Lya PS at different redshifts z and wave¬ 
lengths k. We show k 3 P Lya (k,z)/2n 2 in [nW 2 m - 4 sr -2 ], while k 
is in [/iMpc -1 ]. 


present also for shot noise fluctuations, where fluctuations 
are always suppressed. 

The Fourier transform of eq. ( |65[ ) is 

Sic, Dnt(^k) = L J dn J dz x 

x A(z, z')n u (i/', z')(b(z')) l/ /D(z l )e lk ) S(z = 0,k) = 

= J dz 1 A 1 n v {b')D' ^ (z = 0 , k); (67) 

where in the last equation we avoided to write explicitly the 
dependence in z and z . Any evolution of 2 with x was ne¬ 
glected; this is reasonable because on the scales relevant for 
intensity mapping (10— lOOMpc) the redshift range sampled 
is small. 

The resulting power spectrum is 


k) = 


J dz'A(z, z')(b(z'))„i x 


, / /\ n / ,,sin(fcZ(z,z' 

xn ^’ z)D{z) nd(^ 


/\ "1 2 


al(z = 0). (68) 


As we will see in Sec. |5.4| the IGM emission is very 
smooth on scales smaller than hundreds of Mpc; however it 
is significant if we consider larger scales, comparable to the 
typical redshift distance between the Lya and the Ly/3 line. 
For this reason we neglected shot-noise, because the discrete 
nature of the sources is not important on such large scales. 

Other fluctuations can arise in the Gunn-Peterson op¬ 
tical depth, through the H i density. It affects both the 
absorption probability outside the integral (431 and the 
transmission fractions. These fluctuations are suppressed be¬ 
cause, as discussed in Sec. |4.3[ we can not consider photons 
absorbed in the Lya line since at first order the UV intensity 
field is homogeneous. The bulk of the Lya photons emitted 
in the IGM originates from the absorption of photons emit¬ 
ted between the Lya and the Ly/3 line, therefore the mean 
Lya intensity from the excitations of other Lyman lines is 
dominated by the ISM emission by orders of magnitude. 


5.4 Total power spectrum 

In this Section we will examine the results for the total Lya 
intensity power spectrum predicted by our model. For sim¬ 
plicity we have neglected the correlations among different 
processes and simply summed the individual terms. This 
procedure, although approximate, yields a conservative esti¬ 
mate, as such cross-correlations might increase the power 
spectrum amplitude up to a factor of 2 (since = 

o (a—b)~ \ 

Z a 2 +b 2 >■ 



















12 P. Comaschi & A. Ferrara 


Table [3] shows the numerical values of the tota l Lya PS, 


10 


shows its 


P Lya (k, z), predicted by our model and Fig. 
redshift evolution. As expected the power monotonically de¬ 
creases with redshift, following the analogous mean Lya in¬ 
tensity trend. Since this is the main observable for proposed 
intensity mappers (Silva et al. (20131; Pullen et al. (20141; 


Dore et al. (2014|) with a signal-to-noise too low to perform 


tomographic studies, it is important to understand what 
we can learn from such measurement. The power spectrum 
shape shows two key features: a pronounced change of slope 
at k « 1/iMpc -1 and, for z ^ 6, a bump at k « O.OlftMpc -1 . 
As we will se from Fig. |ll| and Fig. |12| they describe changes 
in the well defined physical properties: shot-noise on small 
scales and IGM emission on large scales. 

In Fig. [IT] we show the contributions to the power spec¬ 
trum from the three processes considered at z = 7. The 
ISM emission dominates the power spectrum by orders of 
magnitude at k ^ 0.1/iMpc” 1 , but on larger scales the IGM 
emission from UV excitations can be significant. The exact 
k in which the PS starts to probe the IGM emission can 
depend on the choices made in our model; however it is im¬ 
portant to understand that IGM emission should be taken 
into account when considering Lya intensity fluctuations on 
large scales and that we can discriminate from an observa¬ 
tion if we are observing the ISM or the IGM emission from 
the slope of the PS. It is also interesting to notice that even 
if in our model ~ 70% of the Lya photons are emitted by the 
IGM, their distribution is very smooth and therefore they 
do not give a significant contribution to the PS in all but 
the largest scales. 

Fig.[l2]shows the PS of the different fluctuations consid¬ 
ered at 2 = 7. At k < 1/iMpc -1 , the PS is always dominated 
by linear fluctuations in the Lya or UV emission from the 
galaxies, therefore even the IGM emission traces the SFR 
distribution. Recombinations in the ionized IGM are always 
subdominant, for the low IGM density. Shot noise is signif¬ 
icant only for k > 1/iMpc^ 1 , but on these scales Lya radia¬ 
tive transfer effects are important, because the Lya emission 
can be spread over an extended halo. Since the shot-noise 
power spectrum is constant over k, any deviation from this 
dependence can be used to constrain these effects. 

Fig. [TT] and [12] refer only to the z = 7 case; at differ¬ 
ent redshifts the qualitative results are similar, while the 
quantitative change can be inferred from Fig. [lO] The only 
difference is that at z < 6 the IGM is not completely opaque 
to Lya photons and therefore the IGM emission from UV 
excitations is less prominent compared to the ISM one even 
at large scales. Finally, in Fig.[l3]we show the redshift evolu¬ 
tion of the PS for k = O.l/iMpc -1 . The dominant feature is 
a smooth decline, i.e. a distinctive, model independent sig¬ 
nature associated with, e.g., reionization is not visible. This 
is expected since Lya intensity mapping is a better tracer 
of SFR history rather than of the IGM thermal/ionization 
state. 

It should be stressed that the power spectrum scales 
as the square of the intensity: for example a factor 3 in 
the Lya escape fraction translates to an order of magnitude 
in the power spectrum. Therefore the quantitative results 
from our model depend strongly on the choices made for the 
poorly constrained parameters that we introduced. However 
the uncertainty on the numerical value of our results does 
not prevent a clear understanding of the physics involved 



Figure 10. Total Lya power spectrum at z = 4, 6, 8, 10, 12 
(from top to bottom). 



Figure 11. Decomposition of the power spectrum in the different 
processes at z = 7. We show the diffuse emission in the IGM from 
UV excitations (red dot-dashed), the emission from the ISM in 
the halos (blue solid); the recombinations in the IGM (yellow 
dotted) and the total PS (black solid). 

in the Lya PS; moreover, the use of an analytical model 
allows a precise control of the dependence of the results on 
free-parameters. 


6 SUMMARY AND CONCLUSIONS 

We have developed an analytical model to estimate the 
Lya emission at z > 4. Our goal is to predict the power 
spectrum of the spatial fluctuations that could be observed 
by an intensity mappig survey. In particular we aim at un¬ 
derstanding what physical processes occurring in the EoR 
can be probed with this technique. Our model uses the latest 
data from the HST legacy fields and the abundance match¬ 
ing technique to associate UV emission and dust properties 
with the halos, computing the ISM and IGM emission con¬ 
sistently. 

The Lya intensity from the diffuse IGM emission 
is 1.3 (2.0) times more intense than the ISM emis¬ 
sion at z = 4(7); both components are fair tracers of 
the star-forming galaxy distribution. However the power 
spectrum is dominated by ISM emission on small scales 
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Figure 12. Contribution to the Lya PS from the different pro¬ 
cesses considered at z = 7. We show total PS (black solid); the 
linear fluctuations (blue solid) and the shot-noise (light-blue solid) 
from ISM emission (eq. ( |49[ ) and ( |55| l) ; from IGM emission from 
UV excitations (red dot-dashed, eq. ( |68| l) ; from recombinations 
in the IGM (dotted lines), from fluctuations in the filling factor 
(thin red, eq. (|64|) ) and in the density (thick yellow, eq. (|60[l). 



Figure 13. Redshift evolution of the Lya PS. We show P L y a 
for k = 0.1/iMpc -1 (solid), k = 0.01/iMpc -1 (dotted) and k = 
1/iMpc -1 (dot-dashed). 

(k > O.OlhMpc^ 1 ) with shot noise being significant 
only above k = 1/iMpc -1 . At very lange scales ( k < 
O.OlhMpc -1 ) diffuse IGM emission becomes important. The 
comoving Lya luminosity density from IGM and galax¬ 
ies, = 8.73(6.51) x 10 40 erg s -1 Mpc -3 and p = 

6.62(3.21) x 10 4 °erg s -1 Mpc -3 at z = 4(7), is consis¬ 
tent with recent SDSS determinations. We predict a power 
k 3 P Lya (k,z)/2-K 2 = 9.76 x 10“ 4 (2.09 x 10 -5 )nW 2 m- 4 sr- 2 
at z = 4(7) for k = 0.1/iMpc -1 . 

The quantitative results from our model depend on the 
choice of the free-parameters. However, at least three points 
appear to solidly emerge from the analysis: (i) Lya inten¬ 
sity mapping is a good probe of LSS and of star formation at 
high redshift; (ii) both the IGM and the ISM produce a sig¬ 
nificant, and potentially detectable, Lya emission, with the 
IGM one being is smoother and becoming important only 
on <1 100 Mpc scales; (iii) the shape of the power spectrum 
is important to constrain the physics of Lya emission and 
to understand the nature of the sources. 


Whether an intensity mapping survey is observation- 
ally feasible is still unclear at this time. Even if the outlook 
for the HI 21 cm line (the first intensity mapping survey 
proposed) seems positive and several instruments are al¬ 
ready producing data or are under construction, other lines 
present completely different challenges (such as line confu¬ 
sion). Therefore further study is necessary, in particular to 
test the extraction of the predicted power spectrum when 
actual foregrounds are included. Moreover, modeling the 
physics involved in the line emission process is an extremely 
challenging task. For this reason, the Lya line emerging from 
galaxies should be carefully calibrated against observations. 
In spite of these difficulties, a successful intensity mapping 
survey would produce such a breakthrough that any effort 
towards this goal is certainly worthwhile. In a forthcoming 
paper we will build upon the present results and assess the 
feasibility of a Lya intensity mapping experiment by pre¬ 
senting possible strategies to control foreground(s) removal. 
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